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Abstract 


The  development  of  practical  numerical  methods  for  simulation  of  partial  differential 
equations  leads  to  questions  of  convergence,  accuracy  (in  time  and  space)  and  efficiency. 
Verification  of  a  computational  algorithm  includes  the  process  of  establishing  a  conver¬ 
gence  theory  for  the  discretized  equations.  It  is  well  known  that  the  long  time  behavior 
of  a  system  may  not  be  captured  even  by  “convergent”  approximating  methods  and 
additional  requirements  must  be  placed  on  the  scheme  to  ensure  the  discretized  equa¬ 
tions  capture  the  correct  asymptotic  behavior.  Even  on  finite  intervals,  there  are  always 
uncertainties  in  the  problem  data  that  can  be  a  source  of  difficulty  for  accurate  simu¬ 
lation  of  nonlinear  problems.  These  uncertainties  lead  to  uncertainty  in  the  computed 
results  and  should  be  considered  as  part  of  the  verification  step.  This  research  gives 
preliminary  results  showing  how  sensitivity  analysis  can  be  used  to  provide  a  practical 
precursor  to  dynamic  transitions  and  quantify  numerical  uncertainty  in  simulations  of 
nonlinear  parabolic  partial  differential  equations. 
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1  Summary  and  Objectives 

Research  Objective:  It  is  well  known  that  the  long  time  behavior  of  a  nonlinear  dynami¬ 
cal  system  may  not  be  captured  even  by  convergent  approximating  methods  and  additional 
requirements  must  be  placed  on  the  scheme  to  ensure  the  discretized  equations  capture  the 
correct  asymptotic  behavior.  This  issue  is  particularly  important  when  one  is  forced  to 
use  numerical  methods  to  evaluate  the  asymptotic  behavior  of  a  closed-loop  control  system 
when  the  mathematical  model  is  defined  by  a  nonlinear  partial  differential  equation  (PDE). 
In  addition,  using  feedback  to  eliminate  or  delay  transition  in  fluid  flows  often  requires  some 
type  of  mechanism  to  predict  that  a  transition  is  about  to  occur.  The  recent  papers  [3],  [4], 
[5],  [21]  and  [22]  provide  considerable  evidence  that,  for  certain  nonlinear  systems  that  occur 
in  fluid  flows,  sensitivity  analysis  can  be  used  to  indicate  a  transition  is  about  to  occur. 
In  [4]  and  [5]  it  was  shown  that  this  information  can  be  used  to  determine  when  to  turn 
on  a  controller  to  prevent  transition.  This  report  illustrates  that  similar  sensitivity  tools 
can  also  be  used  to  provide  insight  into  the  asymptotic  behavior  of  the  closed-loop  system. 
In  particular,  it  is  shown  that  time  varying  sensitivities  can  be  used  to  determine  when  a 
numerical  simulation  is  correctly  predicting  the  longtime  behavior  of  the  response.  In  the 
cases  considered  here,  the  trigger  of  a  transition  can  be  a  known  parameter  (wall  rough¬ 
ness,  etc.)  or  an  un-modeled  uncertainty  in  the  problem  data.  This  includes  uncertainty  in 
parameters,  initial  data,  boundary  conditions  and  forcing  terms.  These  uncertainties  in  the 
problem  data  lead  to  uncertainty  in  the  computed  results  and  should  be  considered  as  part 
of  a  verification  step.  In  addition,  although  we  do  not  address  this  issue  here,  it  has  recently 
been  observed  that  finite  precision  arithmetic  and  sensitivity  to  parameter  uncertainty  can 
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lead  to  orders  of  magnitude  errors  in  simulations  of  simple  nonlinear  convection-diffusion 
equations  (see  [1]  and  [3]).  The  focus  of  this  report  is  on  a  nonlinear  reaction-diffusion 
equations  to  illustrate  the  problem  and  method.  However,  we  first  present  a  simple  ODE 
example  to  illustrate  some  of  the  basic  ideas. 

Collaborators:  Dr.  John  A.  Burns  of  the  Interdisciplinary  Center  for  Applied  Mathematics 
and  the  Department  of  Mathematics  at  Virginia  Tech  served  as  the  main  collaborator  for  this 
work.  The  results  in  this  report,  along  with  the  results  that  are  currently  in  preparation  in 
a  more  extensive  manuscript,  were  the  result  of  an  intense  two  week  period  of  collaboration 
between  the  PI  and  Dr.  Burns  and  a  considerable  amount  of  subsequent  communication. 
In  addition,  the  ideas  for  the  first  example  came  from  the  thesis  problem  and  subsequent 
work  of  Dr.  John  R.  Singler  of  the  Mechanical  Engineering  Department  at  Oregon  State 
University.  His  results  and  observations  for  the  simple  control  problem  example  were  quite 
useful  and  much  appreciated. 


2  Results  of  Funding 

The  following  sections  illustrate  the  key  ideas  of  how  sensitivity  analysis  can  be  used  to 
predict  the  onset  of  transition  in  various  numerical  simulations.  The  first  example  is  that  of 
a  simple  closed-loop  control  system,  and  the  second  section  contains  an  example  of  a  classical 
nonlinear  parabolic  partial  differential  equation.  In  each  case,  by  choosing  to  examine  the 
sensitivity  of  the  state  variable  with  respect  to  a  certain  parameter,  we  are  able  to  produce 
numerical  sensitivity  simulations  which  serve  as  an  indicator  that  a  transition  is  about  to 
take  place  in  the  behavior  of  the  original  state  variable. 


2.1  A  Finite  Dimensional  Example 

We  consider  a  3D  system  that  is  typical  of  those  found  in  the  papers  [4],  [5],  [20],  [21], 
[25],  [26]  and  [27].  However,  we  focus  on  the  role  that  small  constant  disturbances  play  in 
transition  and  illustrate  how  sensitivity  information  can  be  used  to  predict  the  transition 
in  these  cases.  The  system  is  governed  by  three  ordinary  differential  equations  and  has  the 
form 

x(t)  =  A(R)x(t)  +  \\x(t)\\  Sx(t)  +  Bu(t)  +  Dq,  x(0)  =  xq,  (2-1) 

where  A(R)  =  [-^Aq  +  R],  A o  <  0  is  diagonal  and  S  =  —S*  is  skew-adjoint.  In  particular, 
this  three  dimensional  system  is  defined  by 
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where  all  constants  are  positive.  Here,  q  is  considered  a  “small”  constant  (un-modeled) 
disturbance.  For  the  runs  here  we  set  a  =  0.5,  /3  =  0.75,  7  =  1.0,  b\  =  1.0,  62  =  0.5  and 
63  =  0.25.  The  linear  operator  A(R)  is  stable  for  all  R  >  0  and  for  the  no  disturbance 
case  (i.e.  when  q  =  0)  the  dynamical  system  is  also  dissipative.  In  particular,  the  nonlinear 
system  (2.2)-(2.4)  has  a  compact  global  attractor.  The  problem  is  sensitive  to  the  parameter 
q  and  this  sensitivity  plays  an  important  role  in  the  transition  process. 

Let 


A  ax(t,g)  dx(t,  0) 

s(t>  =  =  ~dT 


(2.5) 


denote  the  sensitivity  of  the  solution  x(t)  =  x(t,q)  at  q  =  0.  It  follows  that  the  sensitivity 
s(t)  satisfies  the  linear  differential  equation 


s(t)  =  A(R)s(t)  +  F(x(t))s(t)  +  D,  s(0)  =  0,  (2.6) 


where 


F(x)  = 


||x||  S  +  SxxT  /  ||x||  ,  x/0 
0,  x  =  0 


Consider  the  case  where  x$  =  [9,9, 9]r  x  10-6  and  q  =  5  x  10~8.  Figure  1  below 
contains  the  plots  of  the  norms  of  solution  x(t,q)  and  the  sensitivity  s(t)  (top  plot)  for 
this  system.  Observe  that  the  solution  does  not  “transition”  to  the  (chaotic)  attractor 
until  t  =  175 seconds  .  However,  at  t  =  50 seconds  the  sensitivity  s(t)  satisfies  ||s||  >  103. 
The  vertical  red  line  at  t  =  50 seconds  indicates  that  the  sensitivity  information  provides  a 
precursor  to  the  upcoming  transition  long  before  the  transition  is  observable  in  the  state. 
This  precursor  was  used  by  Singler  to  determine  when  to  switch  on  a  capturing  feedback 
controller  which  is  then  able  to  prevent  the  transition  (see  [20]). 

In  the  next  section  we  use  a  similar  technique  to  investigate  the  numerical  simulation 
of  the  longtime  behavior  of  a  nonlinear  parabolic  PDE.  However,  in  the  PDE  case  the 
“sensitive”  parameter  is  in  the  boundary  condition  which  is  typical  in  parabolic  diffusion- 
convection-reaction  equations  (see  [1],  [3],  [4],  [5],  [13],  [20]  and  [21]). 


2.2  The  Chaffee-Infante  Equation 

We  consider  a  particular  reaction-diffusion  equation  first  studied  by  Chaffee  and  Infante 
in  [9]  and  [10].  This  model  is  a  well  understood  dissipative  dynamical  system  with  a 
global  attractor  consisting  of  a  finite  number  of  fixed  points  and  the  corresponding  unstable 
manifolds  (see  pages  301  -  306  in  [18]  for  details).  In  particular,  we  focus  on  the  partial 
differential  equation 

zt(t,  x)  =  zxx(t,  x)  +  x)  —  [z(t,  x)]2),  0<x<tt,  t>  0  (2-7) 

with  initial  condition 

z(0,x)  =  <f>(x),  (2.8) 

and  Dirichlet  boundary  conditions 

z(t,  0)  =  0  =  z(t,  7r).  (2.9) 

Here  A  >  1,  and  in  this  setting  it  may  be  helpful  to  think  of  (2.7)-(2.9)  as  a  closed-loop 
system  that  we  wish  to  simulate.  It  is  sufficient  to  consider  the  case  where  A  =  4.1  so  that 
the  following  result  holds  (see  page  303  in  [18]). 
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Figure  1:  Norm  of  x(t)  and  s(t ) 


Theorem  2.1.  If  A  =  4.1,  then  the  system  (2.7)-(2.9)  has  five  fixed  points  o(-)  =  0, 

and  cj) 2(-)  inL2(0, 7r).  The  fixed  points  4>o(-)  =  0,  and  4> 2O)  are  unstable 

and  the  attractor  consists  of  the  unstable  manifolds  for  these  fixed  points  along  with  the 
stable  fixed  points  </q  (■)  and 

Figure  2  is  a  schematic  of  the  global  attractor.  However,  for  certain  initial  conditions 
trajectories  are  pushed  rapidly  towards  the  unstable  zero  fixed  point  before  “transitioning” 
to  one  of  the  stable  fixed  points  (•)  or  This  is  similar  to  the  previous  ODE  example 

except  for  the  fact  that  this  system  is  not  chaotic.  However,  if  one  uses  standard  simulation 
schemes  to  investigate  the  dynamic  behavior  of  this  system  it  is  easy  to  obtain  misleading 
results. 

Consider  the  case  where  the  initial  function  is  given  by  fi(x)  =  1.5sin(3x).  Using  the 
MatlabTA/  routine  pdepe  to  simulate  (2.7)-(2.9)  on  the  interval  0  <  t  <  8,  yields  the  solution 
shown  in  Figure  3. 

It  appears  that  by  t  =  2  the  solution  has  “converged”  to  the  zero  steady  state  solution. 
However,  since  the  theorem  above  tells  us  that  this  fixed  point  is  unstable  we  know  that 
this  is  unlikely.  Indeed,  if  one  continues  to  run  the  simulation  to  t  =  16  we  observe  that 
the  solution  actually  “transitions”  to  the  stable  fixed  point  </>]"(•).  This  is  shown  in  Figure 
4  below.  This  is  also  clearly  demonstrated  in  Figure  5  and  Figure  6  which  contain  the  plots 
of  the  L2  norms  of  the  solution  on  the  intervals  [0,8]  and  [0, 16],  respectively. 

In  more  complex  problems  one  does  not  always  have  the  type  of  qualitative  information 
that  is  conveniently  available  for  the  Chaffee-Infante  equation.  As  illustrated  by  this  ex¬ 
ample,  it  is  not  always  clear  when  a  particular  numerical  solution  is  producing  the  proper 
asymptotic  results.  Even  if  the  algorithm  does  eventually  capture  the  correct  limiting 
behavior,  it  is  not  obvious  how  long  one  must  run  the  simulation  to  see  this  result.  There- 
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Figure  2:  Global  Attractor  for  the  Chaffee-Infante  Equation 


Figure  3:  z(t,x)  on  0  <  t  <  8 


fore,  it  is  important  to  devise  numerical  methods  that  can  help  predict  when  a  simulation 
has  “converged”  to  the  correct  asymptotic  behavior.  Equally  important  is  the  ability  of  an 
algorithm  to  generate  and  evaluate  secondary  information  that  might  indicate  when  it  is 
unlikely  that  an  algorithm  has  “converged”  to  the  correct  asymptotic  behavior.  Although 
this  is  a  difficult  problem  for  general  systems,  in  certain  cases  sensitivity  analysis  can  be 
helpful  in  dealing  with  this  issue. 

2.3  Boundary  Sensitivity  for  the  Chaffee-Infante  Equation 

Here  we  consider  the  sensitivity  of  the  Chaffee-Infante  equation  with  respect  to  the  boundary 
condition.  In  particular,  we  replace  the  boundary  condition  (2.9)  with  the  non-homogenous 
Dirichlet  boundary  condition 

z(t,  0)  =  q  =  z(t,  7r),  (2.10) 

where  q  is  a  “small”  number.  It  can  be  shown  that  the  Chaffee-Infante  equation  is  highly 
sensitive  to  changes  in  the  boundary  conditions  and  this  allows  us  to  consider  the  sensitivity 
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Figure  4:  z(t,  x)  on  0  <  t  <  16 


Figure  5:  L 2  norm  of  z(t ,  •)  on  0  <  t  <  8 


variable 


s(t,  x )  = 


a  dz(t,x,q) 


dq 


1 9=0  - 


dz(t,  x ,  0) 
dq 


to  analyze  this  sensitivity  near  <7  =  0.  The  sensitivity  s(t,x)  satisfies  the  linear  boundary 
value  problem 


st(t,  x)  =  sxx(t,  x)  +  X(s(t,  x)  —  2[z(t,  x)]s(t,  x)),  0  <  x  <  ir,  t  >  0,  (2-11) 


with  initial  condition 

s(0,  x)  =  0,  0  <  x  <  7r  (2-12) 

and  boundary  conditions 

s(t,  0)  =  s(t,  7r)  =  1,  t  >  0.  (2-13) 

This  sensitivity  provides  considerable  insight  into  the  long-term  behavior  of  the  solution 
to  the  Chaffee-Infante  equation.  First  note  that  if  the  solution  z(t,x)  —>  4>(x)  and  <j)(x)  is  a 
stable  equilibrium  state,  then  one  would  expect  that  the  sensitivity  s(t,x )  would  approach 


7 


12 


Figure  6:  L 2  norm  of  z(t.  ■)  on  0  <  t  <  16 


a  steady  state  s(x)  satisfying 

0  =  s"(x )  +  A(s(x)  —  2[^(x)]s(x)),  0  <  x  <  7T, 


with  boundary  conditions 


s(0)  =  s(7r)  =  1. 


In  particular,  one  would  have  lim 

t — ^-|-oo 


I  s(t,-) 


\l2 


|s(x)|| 


L2 


c  where  c  is  a  constant. 


This  is  illustrated  in  Figure  7  and  Figure  8  below.  Observe  that  even  though  the  solution 
z(t,x )  appears  to  have  “converged”  to  the  fixed  point  4>o{x)  =  0  by  t  =  2,  and  seemingly 
remains  at  zero  for  2  <  t  <  8  (recall  the  qualitative  information  in  Figure  3  and  Figure 
5),  the  sensitivity  s(t,x )  is  growing  at  an  exponential  rate  on  the  entire  interval  [0,8]  as 
shown  in  Figure  7,  Figure  8  and  is  best  observable  in  Figure  9.  Moreover,  when  the  solution 
transitions  to  the  stable  fixed  point  </>]”(•)  at  t  ~  9  the  sensitivity  is  maximized  and  then 
converges  to  a  (small)  steady  state  as  expected.  If  one  compares  Figure  8  with  Figure  6 
above,  then  it  is  clear  that  this  sensitivity  provides  insight  into  the  transition. 


Figure  7:  s(t,x)  on  0  <  t  <  16 
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Figure  8:  L 2  norm  of  s(t,  •)  on  0  <  t  <  16 


The  most  important  observation  about  the  numerical  results  here  is  that  even  on  the 
“short”  time  interval  0  <  t  <  8,  when  the  numerical  solution  z(t,x )  appears  to  have 
stabilized  at  zero,  the  sensitivity  indicates  otherwise.  In  particular,  in  Figure  9  below  we 
see  the  exponential  growth  of  s(t,x)  on  [0,8]  and  at  t  ~  7  the  norm  of  the  sensitivity 
is  of  the  order  109.  Thus,  even  on  the  short  time  interval  [0,8]  the  sensitivity  provides 
a  clear  indication  that  the  solution  z(t,  x)  has  not  stabilized  at  a  fixed  point  and  that 
it  is  unlikely  that  the  numerical  simulation  is  ” converged”.  As  previously  noted,  this 
insight  can  be  used  to  turn  on  feedback  controllers  to  prevent  transition.  Perhaps  even 
more  importantly,  sensitivity  analysis  of  this  type  can  be  used  to  help  evaluate  numerical 
simulations  in  problems  where  little  is  known  about  the  actual  asymptotic  behavior  of  the 
system. 


Figure  9:  L 2  norm  of  s(t,  •)  on  0  <  t  <  8 
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3  Conclusions  and  Future  Work 


In  this  short  report  we  presented  two  examples  to  illustrate  how  time  varying  sensitivity 
analysis  can  be  used  for  control  and  to  provide  insight  into  the  validation  of  numerical 
simulations  in  nonlinear  systems.  These  ideas  have  also  been  applied  to  a  wide  variety  of 
reaction-convection-diffusion  systems,  and  a  complete  paper  will  appear  in  the  future.  We 
note  that  many  convection-diffusion  problems,  such  as  Burgers’  equation,  show  extreme 
sensitivity  to  boundary  perturbations,  and  sensitivity  analysis  for  these  systems  has  pro¬ 
vided  amazing  insight  into  the  asymptotic  behavior  of  numerical  solutions  (see  [1],  [2],  [3], 
[5],  [6];  [7],  [8],  [12],  [13],  [14],  [15],  [11],  [16],  [17],  [20],  [21]  and  [22]).  Problems  of  this  type 
are  infinitely  sensitive  to  small  parameter  changes  and  can  have  a  dramatic  impact  on  the 
convergence  of  optimal  control  and  design  algorithms. 

We  note  that  although  the  numerical  results  here  provide  considerable  evidence  that 
time  varying  sensitivities  can  play  an  important  role  in  control  design  and  analysis,  consid¬ 
erable  work  needs  to  be  done  to  place  these  ideas  on  a  mathematically  rigorous  foundation. 
We  have  some  theoretical  results  for  parabolic  dissipative  systems  similar  to  the  Chaffee- 
Infante  equations.  However,  much  work  remains  to  be  done.  Finally,  it  is  clear  that  in 
order  to  implement  some  of  these  ideas,  one  needs  to  have  some  indication  of  which  pa¬ 
rameters  (modeled  or  un-modeled)  are  important  to  use  in  the  sensitivity  analysis.  We  are 
currently  looking  into  using  Fisher  information  theory  as  a  mechanism  to  identify  these 
crucial  parameters. 
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